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A Bose-Einstein condensate in a double-well potential features stationary solutions even for 
attractive contact interaction as long as the particle number and therefore the interaction strength do 
not exceed a certain limit. Introducing balanced gain and loss into such a system drastically changes 
the bifurcation scenario at which these states are created. Instead of two tangent bifurcations at 
which the symmetric and antisymmetric states emerge, one tangent bifurcation between two formerly 
independent branches arises [D. Haag et al., Phys. Rev. A 89, 023601 (2014)]. We study this 
transition in detail using a bicomplex formulation of the time-dependent variational principle and 
find that in fact there are three tangent bifurcations for very small gain-loss contributions which 
coalesce in a cusp bifurcation. 
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I. INTRODUCTION 

Bose-Einstein condensates with attractive contact in¬ 
teractions become unstable if the number of particles 
exceeds a certain limit [U3]. In mean-held approximation 
where the condensate is described by the Gross-Pitaevskii 
equation this effect manifests itself in a vanishing ground 
state. Above this limit the mean interaction caused by 
the particles is strong enough to constrict and ultimately 
collapse the condensate. Below this limit the stationary 
solutions of the Gross-Pitaevskii equation are observable 
even though for attractive interactions the ground state 
is not the global minimum of the mean-held energy [4]. 

Characteristic properties of Bose-Einstein condensates 
with attractive interaction could be used for an atomic 
soliton laser [1]. However, this requires the realization 
of a particle how into and out of the condensate. Model¬ 
ing such effects in mean-held approximation is done via 
imaginary potentials, thus rendering the Hamiltonian non- 
Hermitian [5, 6]. Such non-Hermitian systems have been 
widely studied [6-12] and are supported by comparison 
with many-particle calculations [13, 14]. Both in- and 
outcoupling of particles have been experimentally realized 
[15, 16]. 

Despite the non-Hermiticity real eigenvalues and thus 
true stationary states can be found if the Hamiltonian 
is T^T-symmetric [17-19] or pseudo-Hermitian [20-22]. 
This unique property motivated a variety of theoretical 
studies [23-31] and the experimental realization in optical 
wave-guide systems [23, 32-34]. 

In [35] we studied the eigenvalue spectrum and dynam¬ 
ical properties of a Bose-Einstein condensate, described 
by the dimensionless Gross-Pitaevskii equation 

(-A -b Uext + = mV': (1) 


* dennis.dast@itpl.uni-stuttgart.de 


in a three-dimensional T^T-symmetric double-well poten¬ 
tial 


Vext{x) = ^[x'^ + u}'^iy'^ + z'^)] + voe '^'^"-biyxe (2) 


with w = 2, tig = 4, cr = 0.5 and p « 0.12. The symmetric 
real part of the potential is a harmonic trap which is 
separated into two wells by a Gaussian barrier. The 
antisymmetric imaginary part of the potential induces 
particle loss in the left well and particle gain in the right 
well whose strength is given by the gain-loss parameter 7. 

It was found that four stationary states of the real 
double-well potential are created in the two lowest- 
lying tangent bifurcations at strong attractive interac¬ 
tion strengths Na. If gain and loss is introduced into 
the system, instead of two bifurcations there is only one 
bifurcation between two previously independent states, 
while the other two states vanish. However, the under¬ 
lying bifurcation mechanism remained unclear. It is the 
purpose of this paper to clarify this mechanism and study 
the bifurcation scenario in detail classifying it as a cusp 
bifurcation [36]. 

The number of solutions is not conserved at bifurcations 
since the Gross-Pitaevskii equation is non-analytic due 
to its nonlinear part |V'P- This problem is addressed by 
applying an analytic continuation, where the complex 
wave functions are replaced by bicomplex ones. The 
numerical results are then obtained via the time-dependent 
variational principle (TDVP). 

We start with a short introduction to bicomplex num¬ 
bers in Sec. H before discussing their application to the 
Gross-Pitaevskii equation and the TDVP in Sec. HI. A 
detailed analysis of the eigenvalue spectrum and the bi¬ 
furcation scenario is carried out in Sec. IV. Gonclusions 
are drawn in Sec. V. 



2 


II. BICOMPLEX NUMBERS 

Consider a usual complex number z G z = x + 
\y with x,y G M. and the imaginary unit i^ = — 1. A 
bicomplex number is constructed by replacing the real 
and imaginary part of z by complex numbers with an 
additional imaginary unit j which also fulfills = —1, 

z = x+iy = (a;i+ja;j)+i(j/i+j?/j) = zi+izi+jzj+kzk. (3) 

In the last step k = ij with k^ = 1 was introduced. 

Every bicomplex number can be uniquely written [37, 
38] as 

z = z+e+ + z_e_, (4) 

with the components 

z+= (zi + Zk) + i(zi - Zj), (5a) 

z_ = (zi - Zk) + i(zi + Zj), (5b) 

and the idempotent basis 


One can easily confirm that e± fulfill the relations 

e_|_+e_ = 1, e+ —e_ = k, e+e_ =0, = e±. (7) 

Note that an equivalent representation exists in which the 
components z± are complex numbers with the imaginary 
unit j instead of the imaginary unit i. 

The introduction of the idempotent basis significantly 
simplifies arithmetic operations of bicomplex numbers 
since for two bicomplex numbers z,w the relation 

z ow = (z+ o w+)e+ + (z_ o w-)e- (8) 

holds, where o represents addition, subtraction, multi¬ 
plication and division. Thus the bicomplex algebra is 
effectively reduced to a complex algebra within the coeffi¬ 
cients. 


III. NUMERICAL METHOD 

As a next step we apply the analytic continuation 
to the Gross-Pitaevskii equation (1) by allowing bicom¬ 
plex values for the wave function ip and the chemical 
potential y. With the notation introduced in Eq. (4) 
the modulus squared can be written as \p\'^ = pp* = 
p+p'Le+ +p-p*^e-. Here p* denotes a complex conjuga¬ 
tion of p with respect to the complex unit i. Since k = ij 
this complex conjugation changes also the sign of k and 
thus transforms e+ into e_ and vice versa. Due to the 
properties (7) and (8) the ± components are independent 
and the bicomplex Gross-Pitaevskii equation can be split 


into two coupled complex differential equations, which 
evidently do not contain a non-analytical term anymore, 

(-A + Uext + 8'KNap+p*_)p+ = y+p+, (9a) 

(-A - Uext + 8-KNap-p*jP)p- = y-p-. (9b) 

For complex wave functions the relation p^ = p- holds, 
as can be seen in Eqs. (5), and the Eqs. (9) are reduced 
to the Gross-Pitaevskii equation (1). Thus, the coupled 
Eqs. (9) still support all solutions of the Gross-Pitaevskii 
equation (1). 

To solve the Eqs. (9) the TDVP is generalized for 
bicomplex differential equations. Our ansatz consists 
of coupled Gaussian wave packets which yields accurate 
results even for a small number of wave packets [39-41] . 
For the double-well potential studied in this work it was 
shown that two coupled Gaussian wave packets suffice 
and additional wave packets only lead to small corrections 
[35, 42[. For all calculations the following bicomplex 
ansatz is used 

2 

?/>(x, z(t)) = ^ exp(—x^A"x -I- b" • X -I- c"), (10) 

n—1 

with the bicomplex three-dimensional diagonal matrix 
A" = diag(A|(, A" , A” ), the bicomplex three-dimensional 
vector b" = (6”, 0,0)^, and the bicomplex number c". 
The vector z(t) contains all parameters A", b", c" of 
the wave function. The ansatz is chosen such that all 
solutions are invariant under rotations around the x axis, 
i.e. we search only for wave functions that possess the 
symmetry of the potential. The complex ± components 
of the wave function can be written as 

2 

V’±(x,z±(t)) = ^ exp(-x^AJx -I- b!^ • X -I- cj) 

n—1 

(11) 

n—1 

using Eq. (8) since the exponential function is defined as 
a power series. 

The TDVP makes use of the McLachlan variational 
principle [43[ which demands that the variation of the 
functional 

I=\\iy,-npf (12) 

with respect to (f vanishes, then after the variation if = p 
is set. Using the representation with the idempotent basis 
(4) for ip, p and H leads to 

I ={if- - H-p- |i<p+ - 'H+p+)e+ 

+ (i<7’+ - 'H-kV'+liT^- - 'H-p-)e- 
=I+e+ + I-e- (13) 

The two coefficients are complex conjugate, thus, the func¬ 
tional I = 1+ = IZ_ already contains the full information. 
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The variation of the functional I reads, 

51 = {i5^p-\i(p+ - 'H+V'+) + (iv?- - 
= (ii9z_'0-|it^’+ -'H+ij;+)57,*_ 

+ (ii/i_ -'H_'0_|ii9z+'*/'+)<5z+ 

= 0. (14) 

Since the variations 5z^ and Sz'H are independent both 
coefficients have to vanish, which can be written in the 
compact form 

(iipzf — = 0. (15) 

Using the idempotent basis we can reduce the bicomplex 
equation to a pair of coupled complex equations. 

The next step is to insert the ansatz (11) of the wave 
function into Eq. (15). The calculations are in full analogy 
with the TDVP for complex Gaussian wave packets [44] 
and therefore we only list the results. 

The equations of motion for the parameters of the 
Gaussian wave packets read 

iij=4(A^)2-U4, (16a) 

ib^=4A!^b!^+v”±, (16b) 

i4 = 2trA!^-bJ-b!^+z;0±. (16c) 

The coefficients of the effective potential V^, v”j_ and 
no± are obtained by numerically solving 

N 

+•X + <±)gJ) 

n—1 

N 

= + STTiVa^/)* V'±)5±), (17) 

n—1 

with a,f3 = 1,..., d; m = 1,..., TV and fc + Z = 0,1, 2, 
where d is the dimension and N the number of coupled 
Gaussians. For the ansatz (11) used in this work we have 
d = 3 and N = 2. 

Stationary solutions are obtained by varying the pa¬ 
rameters A^, b!^, (c!^ — c^), such that they fulfill the 
equations = 0, b!^ = 0, (c!^—c^) = 0 for n = 1,..., A^. 
Note that the differences {d^ — c^) are used since not all 
parameters c^. are free due to the norm constraint and 
the choice of a global phase. Standard numerical methods 
for complex numbers can be used to achieve this since 
all parameters in Eqs. (16) and (17) are complex. This 
again shows the benefit of the idempotent basis. 

IV. RESULTS 

Before discussing the analytic continuation the real 
spectrum is investigated. In [35] the real eigenvalue spec¬ 
trum of this system has already been studied, however, 
the bifurcation scenario at strong attractive interactions, 
on which we concentrate in this paper, was only discussed 



Figure 1. (Color online) The bifurcation scenario of the VT- 
symmetric stationary states of the three-dimensional double¬ 
well potential at strong attractive interactions. The chemical 
potential of the stationary solutions is shown as a function of 
the interaction strength Na for different values of the gain- 
loss parameter 7 . For 7 = 0 (black dashed line) two tangent 
bifurcations T1 and T2 are present. The states |A) and |D) 
emerge from the bifurcation T1 whereas T2 gives birth to the 
states |B) and |C). For small values of 7 (solid lines) a third 
bifurcation T3 appears at which the two states |B) and |C) 
bifurcate which connects the so far independent branches of 
the states |C) and |D). For higher values of 7 (dotted lines) 
the bifurcations T2 and T3 vanish, as does the intermediate 
state |C), thus, the branches of the states |B) and |D) merge. 


briefly. In particular the bifurcation scenario was only 
shown for the gain-loss parameters 7 = 0 and 7 = 0.001, 
which, as we will see, does not suffice to understand the 
bifurcation process in detail. 

Figure 1 shows the relevant part of the eigenvalue spec¬ 
trum as a function of the nonlinearity parameter Na 
for different gain-loss parameters 7 . All states are VT- 
symmetric, therefore, their chemical potentials /r are real. 
For a vanishing gain-loss contribution, 7 = 0 , the two 
states I A) and |B), which originate from the ground and 
excited state of the double-well potential, emerge from 
independent tangent bifurcations T1 and T2. The bifur¬ 
cation T1 gives birth to the states |A) and |D) whereas 
the states |B) and |G) emerge from the bifurcation T2. 
At 7 = 0.001 the situation has changed fundamentally. 
The two branches |A) and |B) are no longer independent 
but emerge from the same bifurcation Tl. In addition 
the branches |G) and |D) have vanished. 

This behavior can be understood by studying the sys¬ 
tem for very small parameters 7 . For 7 = 0.0004 both tan¬ 
gent bifurcations Tl and T2, and thus, also the branches 
|C) and |D) are still present. However, the two pairs 
emerging from Tl and T2 are now connected by a new 
tangent bifurcation T3 at which the lower lying states 
vanish. If the gain-loss parameter is slightly increased 
the bifurcation T3 is shifted to lower values oi Na and 
approaches T2 (7 = 0.0006). For 7 = 0.0008 both bifur¬ 
cations T2 and T3 have united and vanished. 

To discuss the propagation of the two tangent bifur¬ 
cations in the parameter space in more detail, a phase 
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Figure 2. (Color online) Trajectories of the two tangent bi¬ 
furcations T2 and T3 in the 7 -A^a-parameter space. In the 
shaded area all three states involved in the two bifurcations 
are present while in the other areas only one of these states 
exists. For increasing values of the gain-loss parameter 7 , the 
bifurcation T2 is strongly shifted to lower nonlinearity param¬ 
eters Na. At 7 c ~ 0.000766 the two bifurcations coincide and 
the state |C) vanishes while the states |B) and |D) merge. 


diagram is shown in Fig. 2. The two bifurcations coalesce 
and vanish at critical values of the gain-loss parameter 
7c « 0.000766 and interaction strength {Na)c ~ —1.352. 
Their trajectories have the characteristic form of a cusp 
in the phase diagram. In the area enclosed by the two 
trajectories all three states |B), |C), and |D) involved 
in this bifurcation scenario are present. Outside of this 
area only one state remains. Below the trajectory of T2 
only the state |D) exists and above the trajectory of T3 
only the state |B). Since the areas below T2 and above 
T3 are connected for 7 > 7c, the states |B) and |D) can 
continuously be transferred into each other by following 
a path in the parameter space around the cusp. 

We have seen that the transition from the 7 = 6 case, 
where the states |A) and |B) are completely independent, 
to the 7 = 0.001 case, where both states emerge from 
a common bifurcation, occurs in two steps. The first 
step is the formation of a new tangent bifurcation T3 
between the states |C) and |D) and the second step is 
the coalescence of the two bifurcations. In the following 
we will show that introducing bicomplex numbers to ana¬ 
lytically continue the Gross-Pitaevskii equation allows a 
more detailed discussion of these two steps. 

Bifurcation scenarios are classified by comparison with 
normal forms [36]. The eigenvalue spectrum shows three 
tangent bifurcations which are described by the normal 
form X = — cr. The stationary solutions, x = 0, are 

obviously x = ±i/fT. Two real solutions exist for cr > 0, 
they coalesce at tr = 0, where two complex conjugate 
solutions (cr < 0) emerge. Thus, the total number of 
solutions stays constant if we allow x to be complex. 
Due to the non-analytic modulus squared in the Gross- 
Pitaevskii equation this conservation of solutions is not 
guaranteed although the wave functions can be complex. 

In analogy with the transition from real to complex 
numbers in the normal form of the tangent bifurcation, we 
expect a transition from complex to bicomplex numbers in 
the eigenvalue spectrum of the Gross-Pitaevskii equation. 
By introducing bicomplex numbers the Gross-Pitaevskii 
equation is analytically continued and the two coupled 
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Figure 3. (Color online) Bicomplex bifurcation scenario show¬ 
ing the /ii and component of the chemical potential as a 
function of the nonlinearity parameter Na for different values 
of the gain-loss parameter 7. The i and k components vanish 
for all states discussed, thus fi* — holds and all states are 
VT symmetric. The states with purely real chemical potential, 
i.e. fij = 0, are still present (solid lines). In addition a pair 
of bicomplex states with — fii ± J/rj emerges at the tangent 
bifurcations Tl, T2, and T3 (dotted lines). For 7 > 7c, the 
bifurcations T2 and T3 have vanished leading to a bicomplex 
branch that is completely independent from the remaining 
scenario. Independent of the parameters 7 and Na there are 
four stationary solutions. 


analytic Eqs. (9) are obtained. 

The bicomplex chemical potential fj, of the stationary 
solutions is shown in Fig. 3. For the discussion of the 
results we again use the components of bicomplex numbers 
as introduced in Eq. (3) instead of the representation in 
the idempotent basis since this renders the interpretation 
of the results much clearer. The first thing to note is that 
all states have vanishing /x; and fi]^ components. Therefore 
/X* = /X holds and all states considered are VT symmetric. 
All states discussed so far are still present in the spectrum 
and only have a non-vanishing component (solid lines). 
However, additional states with /xj 0 are found (dotted 
lines). Since the potential in the Gross-Pitaevskii equation 
is symmetric under a complex conjugation with respect 
to j, these states have to appear in pairs /x = /xi +j/ij and 
/X = /xi - j/Xj [ 45 ]. 

First we discuss the case 7 < 7c where all three tangent 
bifurcations are present. At the bifurcation Tl the states 
|A) and |D) coalesce and vanish. On the other side of the 
bifurcation two bicomplex states emerge whose chemical 
potentials are complex conjugate with respect to j. Also 
at the bifurcations T2 and T3 two states of the real 
spectrum vanish and on the other side of the bifurcations 
two bicomplex states emerge. For 7 > 7c the situation 
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Figure 4. (Color online) Comparison of the cusp bifurcation 
occurring in the eigenvalue spectrum (left panels) and in 
the normal form (18) (right panels). For parameter values 
smaller than the critical value (7 < 7c, p < Pc) the tangent 
bifurcations T2 and T3 are separated. They coalesce at the 
critical parameter value (7 = 7 c, p = Pc) and vanish for 7 > 7 c, 
p > Pc. Real branches are drawn as solid lines and (bi)complex 
branches as dotted lines. 


at the bifurcation T1 has not changed, however, the 
bifurcations T2 and T3 have vanished and the states |B) 
and |D) have merged. The bicomplex states that emerged 
from T2 and T3 have also merged and form independent 
branches that do not bifurcate with any real branch. 

At the tangent bifurcations two eigenvectors and the 
corresponding eigenvalues become equal qualifying them 
as exceptional points of second order [46, 47]. The bicom¬ 
plex analysis reveals that at the cusp point, where the 
bifurcations T2 and T3 merge, a total of three states, one 
with real eigenvalue and two with bicomplex eigenvalues, 
coalesce. This is characteristic of an exceptional point of 
third order which has already been observed in spectra 
with similar cusp-like behavior [48, 49[. 

The coalescence of two tangent bifurcations is the char¬ 
acteristic property of a cusp bifurcation [36[ which is 


described by the normal form 

a; = -b px — cr, ( 18 ) 

whose stationary solutions i; = 0 are found by Cardano’s 
method. The spectrum of the normal form has two 
tangent bifurcations which vanish at the critical values 
Pc = o’c = 0. In Fig. 4 the normal form is compared with 
the eigenvalue spectrum of the double-well potential. The 
parameters p and a play the role of the gain-loss param¬ 
eter 7 and the nonlinearity parameter A^a, respectively. 
At the tangent bifurcations in the eigenvalue spectrum 
real values of the chemical potential p turn into bicom¬ 
plex values with a j component, p = pi -I- jpj, whereas 
in the spectrum of the normal form real values become 
ordinary complex numbers. Thus, we have to compare pi 
with Rex and pj with Imx. Again the real branches are 
drawn as solid lines whereas the complex and bicomplex 
branches are drawn as dotted lines. 

In the regime 7 < 7c, p < Pc the two tangent bifurca¬ 
tions are still separated, at 7 = 7 c, p = Pc the tangent 
bifurcations coalesce and finally for 7 > 7c, p > Pc the 
bifurcations have vanished. The normal form captures the 
qualitative behavior of the bifurcation scenario found in 
the eigenvalue spectrum, thus justifying its classification 
as a cusp bifurcation. 

V. CONCLUSION AND OUTLOOK 

This answers the question raised in [35[ regarding the 
bifurcation scenario at strong attractive interactions and, 
thus, completes the discussion of the eigenvalue spectrum 
of the three-dimensional T^T-symmetric double-well po¬ 
tential. Formulating the TDVP with bicomplex numbers 
in the idempotent basis has proved to be useful to analyt¬ 
ically continue the non-analytic Gross-Pitaevskii equation 
and analyze the bifurcation scenarios in the eigenvalue 
spectrum. Using this formalism will help tackling prob¬ 
lems with more complicated PT-symmetric potentials 
and additional interactions such as the dipolar interaction 
which might show an even richer variety of bifurcation 
scenarios. 
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